Flow-Induced Helical Coiling of Semiflexible Polymers in Structured Microchannels 
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The conformations of semiflexible (bio) polymers are studied in flow through geometrically struc- 
tured microchannels. Using mesoscale hydrodynamics simulations, we show that the polymer un- 
dergoes a rod-to-helix transition as it moves from the narrow high-velocity region into the wide 
low-velocity region of the channel. The transient helix formation is the result of a non-equilibrium 
and non-stationary buckling transition of the semiflexible polymer, which is subjected to a compres- 
sive force originating from the fluid-velocity variation in the channel. The helix properties depend 
on the diameter ratio of the channel, the polymer bending rigidity, and the flow strength. 



Buckling is a common phenomena of slender bodies, 
like long filaments and thin sheets, under an external 
load. Examples range from microscopic to macroscopic 
length scales, and include act in filaments, vaulting poles, 
virus capsids, and tectonic plates. Specifically, slender 
rods or semiflexible (bio) polymers exhibit a buckling in- 
stability under compression, when the load exceeds a crit- 
ical value [l|-[3| . In the simplest situation of Euler buck- 
ling of a rod, the symmetry (O2 symmetry along the rod 
axis) is brocken by buckling perpendicular to its axis in 
an arbitrary direction [I,!!!- Since buckling is of such fun- 
damental importance, it has received persistent attention 
over centuries. 

Typically, buckling transitions are considered under 
equilibrium conditions. For macromolecular and bio- 
logical filaments with typical length scales of nano- to 
micrometers, which we are interested in, thermal fluc- 
tuations broaden the buckling transition and generate 
a smooth crossover from the unbuckled to the buck- 
led state [l|-[3|. However, very little attention has been 
a;^d to conformational instabilities far from equilibrium 
1 , which exhibit qualitatively new features compared 
to equilibrium ones. Here, the transport of semiflexi- 
ble polymers in microchannels and capillaries presents a 
new opportunity to study the non-equilibrium behavior 
of such filaments [3-[9| . Vice versa, a detailed understand- 
ing of the dynamical process involved in such a transport 
is of paramount importance in many applications. This 
applies, in particular, to many biologically-relevant poly- 
mers, such as DNA, actin filaments, and microtubules, 
which are semiflexible; an example is DNA sorting in mi- 
crochannels [lol-liij. 

Non-equilibrium instabilities can appear under various 
conditions. Here, we investigate flow fields of spatially 
varying flow strength. Such a situation is easily realized 
in flows through spatially-structured microchannels, e.g.^ 
in a region where the channel width changes from narrow 
to wide (see Fig.[T]). By mesoscale hydrodynamic simula- 
tions of a semiflexible polymer in such a microchannel, we 
observe a buckling of the polymer as it enters the wider 
channel section. Buckling is often the first step in the 
formation of more complex structures [13|, [l^. Indeed, 



FIG. 1: (color online) Polymer conformations in channels with 
the diameter ratios (a) Dw/Dn — 2 and (b) Dw/Dn = 4 for 
the flow strength Pe ^ 550 and L/Lp — 1. See also movies 
SI and S2 in the Supplemental Material [2Q| ]. 



we observe that buckling is followed by a flow-induced 
helical coiling of the polymer. 

Helices are a preferred shape in nature. They are fre- 
quently adopted in biological systems such as proteins, 
they have been shown to be close to optimal packing un- 
der the restriction of their own volume and finite length 
[l3[ , they can arise for entropic reasons in the presence of 
depletant molecules Q, and bacteria swim by rotating 
helical filaments [l5|. Helices are also formed by falling 
ropes and falling fluid filaments lH, [l^ , and transiently 
in homopolymer collapse from an extended conformation 
[Til, or are induced by electrostatic interactions in poly- 
electrolytes [l9|. However, as we will explain, the he- 
lix formation process in microchannel flow is different 
from previously studied cases, because it occurs at small 
Reynolds numbers, in the absence of attractive interac- 
tions, and is mainly due to solvent friction. 

We apply a hybrid simulation method, in which molec- 
ular dynamics simulations for semiflexible polymers are 
combined with a mesoscale simulation method, the multi- 
particle collision dynamics (MFC) approach, for the sol- 
vent [m, [ill. In the MFC algorithm, the fluid is de- 
scribed by a set of N point particles of mass m, with 
velocities determined by a stochastic process. The parti- 
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cle dynamics evolves in two steps. In the streaming step, 
the solvent particles move ballistically for a time interval 
h. In the collision step, particles are sorted into the cells 
of a cubic lattice of lattice constant a and their relative 
velocities, with respect to the center-of-mass velocity of 
each cell, are rotated around a randomly oriented axis 
by an angle a. For every cell, mass and momentum are 
conserved, which leads to the build up of hydrodynamic 
interactions (HI) between the fluid particles; at the same 
time thermal fluctuations are taken into account (2ll-[23| . 

The fluid is confined in a cylindrical channel of diam- 
eter Du,, with periodic constrictions of smaller diameter 
Dn (see Fig.d]). In order to create a smooth laminar flow, 
there is a segment between the narrow and wide parts, 
in which the diameter interpolates linearly between Dn 
and . The length of the constriction is larger than the 
polymer contour length L. No-slip boundary conditions 
are imposed on the channel walls by the bounce-back 
rule and virtual wall particles |22i] . Flow is induced by a 
gravitational force (mg) acting on every fluid particle in 
the direction of the channel axis, which yields parabolic 
flow profiles; in the parts connecting the narrow and wide 
segments, there is also a radial flow component. The dif- 
ference in velocity between the narrow and wide parts of 
the channel is \vn\ — \vw \ = \vn\{l — D'^/D'^) = Av^ where 
Vw and Vn denote the velocities in the wide and narrow 
segments, respectively. 

The polymer is represented by a bead-spring model, in 
which monomers, each of mass M, are linearly con- 
nected by harmonic springs of equilibrium length b 0. 
Excluded volume interactions are taken into account by 
a purely repulsive (truncated and shifted) Lennard- Jones 
potential [22|, |2^ . Additionally, a three-body potential is 
applied to capture bending stiffness, with the bending 
rigidity 0. 

The polymer-fluid interaction is taken into account by 
including the monomers in the collision step [22[ . We em- 
ploy the MFC parameters h = /ksT, a = 130°, 
the mean number of particles per collision cell (N) = 10, 
and the fluid mass density g = {N)m/a^. For the poly- 
mer, we set M = m{N), b = and the Lennard- Jones 
parameters a = b and e/ksT = 1. If not otherwise in- 
dicated, the polymer length is L/b = — 1 = 99. To 
maintain a constant temperature, a local velocity scaling 
algorithm is applied, which yields Maxwell-Boltzmann 
distributed velocities. We characterize the strength of 
the flow by the Peclet number Pe = 7r, where 7 = 
ggDn/{4:r]) is the shear rate at the cylinder wall and r 
the longest relaxation time of a semiflexible polymer [25| . 
Based on a semiflexible polymer model, we find the re- 
laxation time r / Y^7?ia^7^B^ ^ 1.4 x 10^ L/Lp for the 
above parameters, where Lp = K./{kBT) is the persis- 
tence length of the polymer. 

Figure [T] shows snapshots of polymer conformations 
when the polymer leaves the constriction and enters the 
wider part of the channel. Inside the constriction, the 



polymer conformations are rodlike with small fluctua- 
tions. When the polymer is ejected from the constric- 
tion, it enters a region where the fluid velocity decreases. 
Hence, it is subjected to a compressive force as a result 
of the difference in fluid velocity, and the polymer un- 
dergoes conformational changes. With increasing flow 
strength, we find a gradual change in transient polymer 
conformations from (bend) rods to helical shapes. This 
helical conformation is a novel feature for a semiflexible 
polymer, and occurs by spontaneous symmetry break- 
ing. Simulations for various channel diameters D^^ but 
fixed constriction diameter D^, show that the helix di- 
ameter strongly depends on D^, and decreases with in- 
creasing channel diameter ratio D^/Dn, see Fig. [TJ We 
attribute this behavior to the larger compressive force 
due to the larger difference Av in fluid velocities for larger 
Dyj/ Dn- It is evident from the snapshots, in particular in 
Fig. [lb, that the polymers assume helical shapes without 
any direct interaction with the channel walls. Hence, any 
qualitative influence of confinement can be unequivocally 
ruled out as the dominating mechanism. 

In order to understand the underlying mechanisms 
and driving forces, we first study the influence of hydro- 
dynamic interactions. Simulations with turned-off hy- 
drodynamic correlations using the Brownian MFC ap- 
proach where each monomer independently per- 
forms stochastic collisions with uncorrelated fluid parti- 
cles subject to a parabolic flow profile, provide the same 
qualitative behavior. Hence, the rod-to- helix transition is 
not caused by hydrodynamic interactions, although sub- 
stantial quantitative differences in helical properties are 
observed. 

To quantify the helical conformations, we introduce the 
torsional order parameter Ut = Ylf^i~^ cos7i, where the 
angle 7^ is defined by the three subsequent bond vectors 
Ri-i^ Ri, and i^i+i, with 

QQQJ. — (^^-1 ^ ' ^ Ri-\-l) 

\Ri-i X Ri\\Ri X Ri-\-i\ 

Ut is displayed in Fig. [2fa) for the persistence length 
Lp = AL. Ut increases as the polymer forms a helix when 
it is ejected from the constriction, and decreases again as 
the polymer relaxes back to its original straight config- 
uration. Simulations without hydrodynamic interactions 
display a very similar behavior except that the values of 
Ut are significantly larger; this is due to the lack of fluid 
being dragged along by the polymer segments. 

Figure [2fb) presents the maximum values U^ of the 
average order parameter {Ut{z)) for various stiffnesses 
in the range 1/4 < Lp/L < 4. The parameter of more 
flexible polymers increases faster than that of stiffer ones. 
All peak values seem to saturate at high Peclet numbers, 
where the saturation value increases with increasing stiff- 
ness. Within a certain range of Peclet numbers, a univer- 
sal dependence of U^ appears independent of stiffness. 





FIG. 2: (color online) (a) Order parameter Ut as func- 
tion of the axial position of the polymer center of mass 
for systems with (solid lines) and without (dashed lines) 
hydrodynamic interactions. The flow strengths are Pe = 
685,1030,1370,2060 (HI) and Pe = 26,38,51,77 (no HI), 
the ratio L/Lp = 1/4, and D^/Dn = 2. (b) Peak values 
Ut of the order parameter Ut for the persistence lengths 
L/Lp = 4.0 (□), 2.0 (O), 1.0 (A), 0.5 (<), and 0.25 (v). 
The results are for the MPC fluid with Du^/Dn = 2. 



In order to elucidate the factors which cause helical 
coiling of a polymer under a non-homogeneous compres- 
sive forces, we consider a simpler system, where a poly- 
mer in a viscous fluid is driven by a rigid wall (penetra- 
ble to the fluid), which moves with the constant velocity 
v^. Thermal fluctuations and hydrodynamics interac- 
tions are omitted for simplicity. The planar wall is par- 
allel to the x?/-plane and perpendicular to the polymer, 
which is oriented initially along the z-direction with one 
end in contact with the wall. The polymer dynamics is 
governed by 

Mr, = -Tr, + F^^Giz^ - z,) ^ ^ (2) 

where is the position of monomer the position 

of the wall, 0(z) the Heaviside step function, Fl and F^ 
are the bond and bending forces, respectively, and T is 
the monomer friction coefficient. F^ = Fe^, where 
is the unit vector along the 2;-axis, describes the force 
of the wall on monomer i, when it is in contact with the 



FIG. 3: (color online) Typical configurations of buckled semi- 
flexible polymers with Lp/L — 1.0. The rigid substrate is 
indicated by the gray area. The polymer length is Nm = 100 
and a the force is F / {UbT / a) — 100. The viscous drag coeffi- 
cient is (a) V I ^ mkBT jh^ = 1.0 and (b) 50.0. The sequence 
of images (c) - (d) illustrating the helical coiling of the semi- 
flexible polymer during the initial phase. See also movies S3 
and S4 in the Supplemental Material [ 20| 



wall. Equation [2] is solved numerically for a polymer with 
Lp/L = 1. 

For small F, less than a threshold frictional dissi- 
pation is weak and the polymer moves with the wall with- 
out any significant conformational changes. For fixed v^, 
the behavior changes qualitatively when F exceeds 
because the compressive force due to polymer friction 
now leads to a buckling transition, similar to classical 
Euler buckling. As the wall moves further, the part of 
polymer in contact with the wall aligns with the wall 
and moves along the wall, whereas the part further away 
from the wall moves in the same direction, but a smaller 
velocity (cf. Fig. [3ti). Ultimately, the entire polymer is 
stretched out parallel to the wall. With further increas- 

(2) 

ing F, beyond a threshold Fc , we observe a distinctly 
different behavior, as the polymer now exhibits a tran- 
sition to a coiled state (cf. Fig. [Sb). As displayed in 
Fig. m the bending energy Ub of the polymer exhibits 
two different regimes as function of time for such systems. 
Initially, Ub increases by local Euler-like buckling of the 
polymer over a length scale much smaller than its length 
(cf. Figs. [31^ c), a)). When the bending energy exceeds 
a certain value, the polymer undergoes a conformational 
change to a coiled state near the wall (Fig. [3l^d),(e)), and 
Ub drops significantly. For later times, the polymer con- 
tinues coiling until the entire chain is coiled up at the 
wall. 

The dependence of the peak value of the bending en- 
ergy on the wall velocity and the friction coefficient F 
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FIG. 4: (color online) (a) Time dependent bending energy 
Ub of a semiflexible polymer with L/Lp = 1.0 pushed for- 
ward by a rigid wall moving with velocity Vw — 0.1 y^sT/m. 
Curves from left to right correspond to decreasing drag co- 
efficients V/^mkBT/h^ = 100,50,20,10,5. (b) Dependence 
of the maximum bending energy JJJ^"'^ (top), the minimum 
energy after the drop JJJ^'^^ (middle), and the minimum bulk 
energy UJ^'^^ (bottom) on friction. The slopes of the straight 
lines are 1/3 (top, bottom) and 1/2 (middle). 



can be understood from some scaling arguments. For 
F < V^c^ , the motion of the polymer through the viscous 
medium generates a inhomogeneous (negative) tension 
(j{s) = Vvw{L — s), there s is the arc length measured 
from the polymer end at the wall. For a homogeneous 
tension ctq buckling theory predicts an instability with 
a fastest growing mode with wave vector q ^ y^ao/tz^ 
hence buckling occurs when q > 27r/L. If we approxi- 
mate the tension by its average <Jo = we obtain 



.(1) 



(2) 

For friction coefficients F > Fc , initial buckling occurs 
at the near- wall end of the polymer (Fig. El^l)). The 
bending energy of the buckled part of contour length Lc 
has the local curvature 1/R{s) and, hence, the bending 
energy is Ub ^ k, J^"" R~'^ds. Since there is only local 
bending of the polymer over a length scale Lc ^ we 
find Ub ~ n/Rhy replacing R{s) by a characteristic value 
R. The resulting bending force = —n/R^ is balanced 



by the frictional force over the same length scale, i.e., 
Ff = Tvy^R^ which implies a curvature radius 



and a bending energy 



(3) 



(4) 



During the transition from local buckling to circu- 
lar coiling, the curvature of the buckled part decreases, 
which explains the drop in bending energy (Fig. [3^2), (3)). 
At later times, during the stationary coiling process, the 
dynamics is determined by two contributions: the bend- 
ing energy of a part Lc not in contact with the wall 
[755 (bulk part) and a part of circular shape at the wall. 
For the bulk part of the bending energy, the same argu- 
ment applies as in the derivation of Eqs. (j3]), (j4]). Hence, 
we expect to find the same scaling relation for bending 
in the bulk. 

Figure H^b) shows the scaling behavior with respect 
to friction of the maximum bending energy U^^^ of 
Fig. m^a). This energy nicely follows the scaling predic- 
tion (|4]). To characterize the the quasi-stationary behav- 
ior of the bending energy for long times, we determine the 
minimum energy /7^*^ of the bulk part of the bending 
energy after the steep drop of The numerical simula- 
tions show that the bulk bending energy dominates over 
the wall part for small friction coefficients. At the larger 
F the two contributions become comparable. As shown 
in Fig. HKb), UJ^'''^ also follows the scaling relation (g]), 
whereas the total minimum energy exhibits the stronger 
dependence F^/^. 

The buckling and helical coiling of a polymer under 
flow in a structured microchannel and a polymer pushed 
forward by a wall in a viscous medium is governed by 
very similar mechanisms. In contrast to the moving- wall 
model, the polymer in the structured microchannel does 
not exhibit circular collapsed conformations. Indeed, the 
compressive force due to the fluid velocity difference is 
less restrictive, because the coiled part of the polymer 
can relax more easily. The observed behavior differs 
from classical Euler buckling in several respects. First, 
the buckling force appears due to a non-equilibrium and 
non-stationary dynamical process. Second, buckling is 
initiated locally. This is a consequence of the frictional 
force, which increases with the length of the considered 
polymer part. Third, a transition to helical coiling occurs 
when the local tilt angle exceeds a certain value. 

Thus, helical coiling is expected to be a generic phe- 
nomenon of polymer motion in non-homogeneous vis- 
cous environments. Act in filaments [H in structures mi- 
crochannels seem to be an optimal system to observe 
the helical-coiling transition experimentally. 
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